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ABSTRACT 

Context. Extragalactic jets are inferred to harbor dynamically important, organized magnetic fields which presumably aid in the colli- 
mation of the relativistic jet flows. We here explore by means of grid-adaptive, high resolution numerical simulations the morphology 
of AGN jets pervaded by helical field and flow topologies. We concentrate on morphological features of the bow shock and the jet 
beam behind the Mach disk, for various jet Lorentz factors and magnetic field helicities. 

Aims. We investigate the influence of helical magnetic fields on jet beam propagation in overdense external medium. We adopt a 
special relativistic magnetohydrodynamic (MHD) viewpoint on the shock-dominated AGN jet evolution. Due to the Adaptive Mesh 
Refinement (AMR), we can concentrate on the long term evolution of kinetic energy dominated jets, with beam-averaged Lorentz 
factor r =a 7, as they penetrate into denser clouds. These jets have near-equipartition magnetic fields (with the thermal energy), and 
radially varying T(R) profiles within the jet radius R < Rj maximally reaching T ~ 22. 

Methods. We use the AMRVAC code, employing a novel hybrid block-based AMR strategy, to compute ideal plasma dynamics in 
special relativity. We combine this with a robust second-order shock-capturing scheme and a diffusive approach for controlling mag- 
netic monopole errors. 

Results. We find that the propagation speed of the bow shock systematically exceeds the value expected from estimates using beam- 
average parameters, in accord with the centrally peaked T(R) variation. The helicity of the beam magnetic field is effectively trans- 
ported down the beam, with compression zones in between diagonal internal cross-shocks showing stronger toroidal field regions. 
In comparison with equivalent low-relativistic jets (r ^ 1.15) which get surrounded by cocoons with vortical backflows filled by 
mainly toroidal field, the high speed jets demonstrate only localized, strong toroidal field zones within the backflow vortical struc- 
tures. The latter are ring-like due to our axisymmetry assumption and may further cascade to smallscale in 3D. We find evidence for 
a more poloidal, straight field layer, compressed between jet beam and backflows. This layer decreases the destabilizing influence of 
the backflow on the jet beam. In all cases, the jet beam contains rich cross-shock patterns, across which part of the kinetic energy 
gets transferred. For the high speed reference jet considered here, significant jet deceleration only occurs beyond distances exceeding 
O(l00Rj), as the axial flow can reaccelerate downstream to the internal cross-shocks. This reacceleration is magnetically aided, due 
to field compression across the internal shocks which pinch the flow. 

Key words. ISM: jets and outflows - Galaxies: jets - methods: numerical, relativity 



1. Motivation 

Relativistic jets represent extremely energetic phenomena in as- 
trophysics. They are associated with (1) a variety of compact 
objects, (2) with Active Galactic Nu clei (AGN) carrying en- 
ergy fluxes of 10 45 - 10 48 ergs/s ( Celotti et al.|l997 Tavecchio 
|et al.||2TjQ0"| ), or (3) with micro-quasar systems (energy flux of 
10-"ergs/s). These high energies are somehow extracted from 
the inner part of the system and this energy is transported over 
long distances by means of a detectable collimated jet. A large 
amount of this energy gets deposited by the jet in the surround- 
ing medium, as only a small fraction of the jet energy is dissi- 
pated in the innermost region (Sambru na et al.|2006a[ ). In many 



AGNs, the kinetic energy flux is comparable to the bolometric 
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radiative luminosity of the central part (Rawlings & Saunders 
[T991] |Xu7Livio, & Baump999l [Sambruna et al.||2006a) . This 
implies that a study of the interaction of the magnetized jet with 
external medium could help provide model constraints for jet 
formation and collimation. The jet is assumed to be powered 
by a spinning black hole (Blandford & Zn ajek|1977| |Begelman| 



et al. 1984|i, and/or b y the surrounding accretion disk corona 
(Mill er & Stone|2000[ ). The disk-jet is launched by general rel- 
ativistic magnetohydrodynamic (GRMHD) mechanisms, and is 
accelerated to high Lorentz factor. The observations show that 
AGN jets propagate in the parsec scale with a Lorentz factor 
F ~ 10-30 ( |Kellermann et al.|2004) , and in some AGN types 
(blazar, QSO), the jets are relativistic even at kpc scale with 
F ~ 5-30 ( |Tavecchio et al.|2004[ ). In such objects, the radio toX- 
ray observations cannot be explained by the presence of a unique 
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synchrotron component. Scenarios involving synchrotron and 
inverse Compton scattering of synchrotron photons are able to 
fit observations, under the assumption that the jet decelerates 
from the subpc scale to the larger scale ([Georgeanopou los &| 



Kazanas 2003). This can be achieved if an outer slower layer 
surrounds the central relativistic jet "spine" (|Ghisellini et al 



2005 ), since continuous electron acceleration acts at the bound- 



ary layer between jet sheath and "spine" ( Stawarz & Ostrowski 
2002) l. This kind of jet structure is supported by direct radio ob 



servations ( Giroletti et al.|2004 l and leads to radial jet structure 



where the Lorentz factor decreases radially. Such a stratification 
of the jet where the Lorentz factor increases towards the axis, has 
been suggested from observations, both in microquasars ([Meier 
[2003) , and in AGN Pester et al.|2006||Dulwich et al.|2007) . Such 
transverse jet stratification is also supported by models of jet 
launch scenarios ( |Koide et~ 2001 ; McKinney 2006; Meli ani et| 
|al.|2006] >. However, most numerical investigations to date have 
ignored radial stratification of the jet velocity, and have made 
simplifying assumptions on the magnetic topology. In this paper, 
we present numerical simulations of the propagation of magne- 
tized, relativistic, radially stratified jets. These are important to 
understand the impact of more realistic jet topologies on the sur- 
rounding medium, and how this in turn helps to deduce proper- 
ties of jet formation mechanisms. We look in particular to the 
effects of helical magnetic fields on the relativistic jet propaga- 
tion in external media. 

The magnetic field seems to play a significant role in the jet 
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2006). There is observational evidence for an intrinsic magnetic 



field in relativistic jets. In BL-Lac objects, VLBI observations of 
the polarisation of the synchrotron emission, and especially the 
rotation measure, have shown that the magnetic field is system- 



atically tilted from the jet axis ( |Asada et al.|2002||Gabuzda et al. 



|2004| l, which indicates the presence of helical magnetic fields. 

Significant progress has been made regarding the numerical 
modeling of relativistic jet propagation, especially for relativistic 
hydrodynamic models (Hard ee"et al.|2005||Perucho et al.|2 006). 
Relativistic MHD studies of jet propagation have been under- 
taken by several authors, mostly restricted to axisymmetric (or 
2.5D) configurations. Recently, 3D simulations of magnetized 
'spine-sheath' jets have emerged as well ( jMizuno et al.||2007) , 
for purely axial magnetic field configurations at modest Lorentz 
factors of T = 2.5. Combined with linear stability results for 
relativistic 'top-hat' profiles (i.e. uniform jet bounded by a uni- 
form sheath), these results confirm the possibility of jet stabi- 
lization by invoking radial structure. We will restrict our study 
to axisymmetric configurations (as even our grid-adaptive simu- 
lations are still fairly computationally expensive for long-term 
evolutions), but emphasize high Lorentz factor flow regimes 
in helical field topologies. Original RMHD simulations were 
presented by |van Putten] ( |1996| l, where low Lorentz factors in 
toroidally magnetized jets were simulated. A significant step for- 
ward was presented by Komi ssarov| ( |1999| , with studies of light, 
Lorentz factor r = 10 jets pervaded by purely toroidal mag- 
netic fields. By confronting a Poynting flux dominated, with a 
kinetic energy dominated jet, clear morphological differences 
were identified and explained: extensive cocoons form in ki- 
netic energy dominated cases, while 'nose cones' develop with 
high stand-off distances between Mach disk and bow shock 
for the Poynting dominated jets. The latter correlate to highly 
magnetized, shocked plasma, self-collimated by strong magnetic 
pinching. These two distinctive cases have been used subse- 



quently by various authors, also as benchmarks for the rapidly 
growing code pare capable of performing special relativistic 
MHD computations (Keppens & Meliani 2007, Mignone et al. 
|2005[ ), while a first more comprehensive study, including also 
purely poloidally magnetized jets was presented by |Leismann et 
al. (2005 ). The latter confirmed the purely toroidal field cases, 
but also demonstrated that purely poloidal field cases did not de- 
velop pronounced nose cones. The latter were found less suscep- 
tible to internal shock deformations, due to the stabilizing mag- 
netic tension. Hence, fairly smooth cocoons and stable beams re- 
sulted. The authors speculated how the resulting brightness dif- 
ferences would allow to distinguish inherent magnetic topolo- 
gies. We here revisit this suggestion, by making a more con- 
centrated effort on kinetic energy dominated jets alone, with a 
more gradual change from toroidal to near poloidal fields. The 
helical field configurations presented indeed confirm the ear- 
lier trends, but our grid-adaptive models do show rich shock- 
structured beams also for more poloidal field regimes. The code 
used here is the AMRVAC code (|Keppens et al.|2003 Meliani et 



al.|2007 van der Hoist & Keppens 2007 1, which has also been 
tested on a variety of stringent ID and multi-D relativistic MHD 
problems ( van der Hoist et al7||2008) l. In section [2j we list the 
governing equations, Sect. [3J provides details on the numerical 
strategy, while Sect.|4]discusses the main results. 



2. Magnetohydrodynamics in special relativity 

In special relativistic theory, where material particles move 
through four-dimensional spacetime with speeds strictly less 
than the speed of light c, the governing conservation laws ex- 
press particle number conservation, and energy-momentum con- 
servation expressable as the vanishing divergence of a stress- 
energy tensor. This latter four-tensor includes (perfect) gas as 
well as electromagnetic contributions, and in general, the full 
dynamics governed by gas and electromagnetic field variations 
needs to solve the full set of Maxwell equations as well. When 
we adopt the ideal MHD approximation, in which the electric 
field in the comoving frame vanishes identically, a numerically 
convenient set of conservation laws results when we choose a 
'lab' Lorentzian reference frame with four-coordinates (ct,x) T , 
where x = (xu X2, x%) are the three spatial orthogonal coordinate 
axes. The spacetime metric for special relativistic applications 
is the usual Minkowski metric g a & = diag(-l, 1,1,1), where the 
greek indices take values from 0, 1, 2, 3. If we indicate the proper 
density as experienced in the local rest frame with p = mo«o, 
where niQ denotes the particle rest mass and «o the proper num- 
ber density, Lorentz contraction results in a lab frame number 
density F«o, with the Lorentz factor T = 1/ ~\/(l-v 2 /c 2 ). Particle 
number conservation as written for the fixed reference frame is 
then written using the variable D = Tp as 



— +V-(Dv) = 0. 

at 



(1) 



This includes the three-velocity v, whose magnitude v < c. The 
temporal component of the divergence of the stress-energy can 
be written in the same laboratory frame as 



d_ 
Jt 



(t + Dc 2 ) + V • S = . 



(2) 



The energy flux S expression will be given below, and the total 
energy in the lab frame is then split off in a rest mass contribu- 
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tion Dc 2 , and the partial energy r made up of gas, magnetic and 
electric field energy densities from 



In Eq. B and E indicate the magnitude of the usual three- 
vector magnetic B and electric E fields, respectively. When we 
adopt a simplifying polytropic equation of state where the spe- 
cific internal energy of the gas is directly related to the proper 
(rest frame) density and pressure, i.e. p/(y - l)p with polytropic 
index y, the gas contribution is 



T S =P 



c 2 + 



yp 



(7 -Dp 



r 2 - P - Dc 2 



(4) 



The term between brackets represents the relativistic specific en- 
thalpy containing a rest mass contribution. The energy flux in 
Eq. Q can also be split in a gas and electromagnetic contribu- 
tion as in S = S g + S em , where the latter are then quantified from 



S g = p\c 2 + 

*S - 



yp 



" r \~ (r-i)p 

ExB 



Ho 



(5) 



Obviously, S em is the Poynting flux. With this notation, the spa- 
tial part of the stress-energy divergence, using the ideal MHD 
approximation, can be written in the lab frame as 



as 

~dt 



+ V- Sv + p tot c 2 I- 



c 2 BB (v B) 



j"0 



r 2 



j"0 



vB UO. 



(6) 



The I denotes the three by three identity matrix, and the total 
pressure p lot is computed from 



Plot = p + Pa 



1 /B B (v-By 



2p \ r 2 



(7) 



The system is then closed with the homogeneous Maxwell equa- 
tions 



V B = 0, 
<9B 



dt 



V x (v x B) . 



(8) 



The latter is completely identical to the non-relativistic induc- 
tion equation in ideal MHD formulations, since vanishing elec- 
tric fields in the comoving frame implies similarly E = -v x B. 
This allows to write the electric energy density in Eq. Oil as 
€qE 2 /2 = e Q [B 2 v 2 - (v • B) 2 ]/2, and the Poynting flux as S7 m = 
[B 2 v-(vB)B]M>. 



where (1q = 1 = £q. This scaling will be exploited from here on- 
wards. In terms of the conserved variables [D,S,t, B] the usual 
non-relativistic ideal MHD equations in terms of the conserved 
quantities [p,m = pv,e = pv 2 /2+p/(y-l)+B 2 /2,B] are directly 
obtained from the limit F — > 1 . In our conservation-law oriented 
integration strategy, the induction equation dS} is written as 



on 

~dt 



+ V ■ (Bv - vB) = % V(V ■ B) . 



(10) 



The coefficient rja is chosen to correspond to the maximal al- 
lowed diffusion coefficient which still complies with an unmodi- 
fied Courant-Friedrichs-Lewy constrained time step At. This acts 
to diffuse potential numerical monopole errors at their maximal 
rate, and has been used in various non-relativistic MHD appli- 
cations. We refer to Keppe ns" et al.| ( |2003| l for a comparitive 
study between this and other popular source term strategies for 
V • B treatments in an AMR framework, and to | van der Hoist &] 
|Keppens| ( |2007| l for its use in AMR in combination with curvi- 
linear coordinates. 

When we further introduce the auxiliary variable ^ from 



? = r> + 



yp 
(r-D 



(ii) 



the Lorentz factor in essence depends on F(S,B;£), as one 
can write v = (S + £ _1 (S • B)B) /(£ + B 2 ). The defining rela- 
tion Eq.Q then becomes 



y 



r 2 



+ B £ - 



l 



r 2 



(S • B) 2 



-t-D = 0. (12) 



A root-finding algorithm (Newton-Raphson) is implemented to 
compute £ from this expression. This then suffices to make 
the needed conversions from conservative to primitive variables 

[p,y,P,B]. 

The shock-capturing conservative discretization (with the 
monopole source term strictly speaking destroying perfect con- 
servation for B only) is then a TVDLF ( |T6th & Ods1rHfT996] l 
type method, in which only the maximal physical propagation 
speed needs to be computed. This is done using a (slightly mod- 
ified from Numerical Recipes) Laguerre's method to compute 
both speed pairs pertaining to slow and fast magneto-acoustic 
waves as the roots of a quartic polynomial (see e.g. |del Zanna et 
[in ([2003)) . Since all 4 roots A must lie in the interval ] - 1, +1[ 
but can become notoriously close to each other and unity, it helps 
to transform (Berg mans et al.|2005[ ) to a variable p = 1/(1 — A) 
with well-seperated roots on ]0.5, +oo[. The TVDLF scheme is 
used with a Hancock predictor, and is second order accurate for 
smooth solutions. This implies limited linear constructions to 
obtain cell edge from cell center quantities, and in this process 
we employ the set [p, Tv, p, B]. 



3. Numerical approach 

3. 1 . Algorithmic details 

For the numerical solution of the relativistic MHD equations, 
we actually combine Eq. (|2| and the particle conservation law 
Eq. ([T]i to obtain the following conservation law 

|+V.((t + Ao 0v-^bU. (9) 
dt \ po ) 

It is then also necessary for the numerical approach to exploit 
a scaling where c — 1 and thus setting electromagnetic units 



3.2. AMR strategy and numerical setup 

All jet simulations are done a domain size [R, Z] e [0, 40] x 
[0,200], and allow for six grid levels (including the base level), 
achieving an effective resolution of 3200 x 8000. The jet internal 
structure will be specified in detail below, and it initially occu- 
pies only the region [Q,Rj] x [0, Zj]. The exterior region will al- 
ways represent a higher density medium in our simulations and 
is meant to mimick a denser 'cloud', leading to jet deceleration 
by entrainment and dissipation, at least in kpc scale. Our nor- 
malization always sets the jet radius Rj - 1.5 and height Zj = 3. 
We will typically compute till times beyond t = 210, and due to 
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the adopted scaling where the light speed c — 1 , we thus follow 
jet propagation for at least 70 light crossing times of the jet beam 
diameter 2R y To convert from dimensionless, computed values 



to physical quantities, one may adopt (Harris & Krawczynski 
|2006| l a typical radius for an AGN jet of 0.05 pc, a cloud num- 
ber density of 10 cm -3 , and the light speed. Our AMR scheme 
exploits a Richardson type error estimator to dynamically cre- 
ate or remove finer level grids where needed. In addition to this 
estimator, we always enforce the highest grid level around this 
inlet corner, such that we have 120 grid cells through the jet ra- 
dius Rj and the same amount of cells through the initial jet height 
Zj. In the Richardson process itself, two low -order f" +1 solutions 
with coarsened grid spacing 2Ax are constructed from the known 
solutions at resolution Ax at times f and f -1 , by reversing the 
order of the time integration and coarsening operations. When 
a weighted average of selected components exceeds a tolerance 
parameter e to i = 0.005, new grids are created. In all runs, these 
selected components include the lab 'density' D, partial energy 
density t, and magnetic field component B v , with weight ratios 
2:1:1. 

The boundary conditions for the simulations enforce the 
primitive variable profiles as discussed in the next section within 
R < Rj at the lower boundary. The remainder of this bot- 
tom Z = boundary is treated as a symmetry boundary for 
D,S r,S v ,t and Bz while we ensure the vanishing of Sz, Br, B^. 
This acts as a kind of reflecting underlying 'disk' configura- 
tion from which potential backflows deflect. The lateral and top 
boundaries are open, while the symmetry axis R = enforces 
the usual (a)symmetry conditions. 



4. AGN jet computations 

4. 1 . Jet inlet conditions 

In all 8 simulations discussed below, we fix the polytropic index 
7 = 5/3. In fact, all models investigated in this paper represent 
underdense jets, and have typically maximum Lorentz factor at 
the axis of 22. Then, we can expect from analogous ID Riemann 
problems that the forward shock and reverse shock are both near 
Newtonian, so that it is quite adequate to use this Newtonian 
polytropic index value. The magnetic configuration in the lab 
frame is initially given by 



B 7 



2B 



B c 



^(f) 3 tanh(|) 4 tanh(|) 2 



§■ cosh 
Bo 



(I) 



cosh(jf) cosh(|) 
= fBitanh(f) in [0,/?j] x [0,Zj], 



elsewhere. 



(13) 



This solenoidal field is fully helical internal to the jet whenever 
B\ + 0, while it is purely poloidal elsewhere. The parameter "a" 
for the azimuthal field variation is held fixed at a — 5, and in all 
but our purely toroidal field simulation (where the exterior is un- 
magnetized), the magnetic field strength of the cloud B c = 0.01, 
corresponding to a weak uniform background magnetic field. We 
explore the relative importance of the jet internal toroidal versus 
poloidal field components by varying the parameters Bo and B\. 
The proper densities of the jet and cloud are typically fixed at 



Pj = 100 and p c = 1000, except for one run where the cloud den- 
sity p c is decreased twofold. Hence, we restrict the discussion to 
the propagation of relativistic underdense jets, since they yield 
deceleration of the jet and formation of the typically complex 
cocoons. Moreover, we explore only relativistic jets dominated 
by kinetic energy, which are then characterised by a strong shock 
(Appl & Camenzind 19 88]), which i n turn are efficient in parti- 
cle acceleration ( |Begelman & Kirk||l990| l. The initial pressure 
distribution is computed from 



P = P) + \bI~\ (Bl(R,Z) + B 2 Z (R,Z)) , 



(14) 



where the reference value for the parameter pj = 1 (representing 
a cold jet). Finally, the initial flow field v is quantified by the 
following prescription 



B w 



v z - a 
B, 

VR = 0. 



y/p\R/a)' 



(15) 



It vanishes outside the jet region [0, Rj] x [0,Zj], and includes 
jet rotation whenever B\ + 0. Note that the requirements p > 
and v < 1 restrict the possible choices for the parameters p s and 
a, once values for Bo, B\ (B c ) and pj are adopted. This prescrip- 
tion of the primitive variables p, v, p and B is inspired by pre- 
vious non-relativistic ideal MHD models as exploited in Casse 
& Marcowith (2005j), and is such that it then ensures the radial 
force balance along the lower boundary Z = 0. We quantify in 
what follows how these parameters translate into various dynam- 
ically important dimensionless parameters like Lorentz factor, 
Mach number, Alfven Mach number, magnetization parameter, 
as averaged over the jet radius. In our implementation, we found 
it convenient to prescribe these jet profiles at the inlet bound- 
ary Z = directly from the above primitive variable expres- 
sions within those bottom ghost cells with cell centers Z < and 
R<Rj. 

In Table [T] we list the input model parameters for the eight 
cases studied. The first three runs only differ in the value of 
a, and mimick the effect of going from nearly non-relativistic 
speeds (NR), to mildly relativistic (MR), to fully relativistic 
flows in our reference case (Refl). The next three models (Pol, 
Tor, Ref2) explore variations in magnetic field topology (from 
almost purely poloidal to fully toroidal). The final two cases 
look at modest changes (factor 2) in external density (Ref3) 
or overall pressure (Ref4). Note that for all the cases consid- 
ered, the entire jet rotation profile enforced at the inlet Z = 
is within the light cylinder, i.e. within the jet radius R < Rj, we 
have | - vzB^/Bz | increasing from zero to 0.63 for most 
of the models. In Table [2] we list the corresponding dimension- 
less parameters which characterize the internal jet properties. We 
list the Lorentz factor F, the Mach number M - vz/c s with 
the sound speed given by c s — -Jy(y — l)/>/((y — l)p + yp), 
and the relativistic proper Mach number M - MT/r s where 
F? = 1/ V(l ~ c s) is the Lorentz factor computed from the sound 
speed. We further give the Alfven Mach number M a = vz/c a 
where the Alfven speed is found from ( Lichnerowicz| 1967 1 



' cr + 1 + 



yp ' 

(y-Dp 
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The latter contains the ratio of magnetic to rest mass energy den- 
sity 



cr = 



2p, 



mag 



which together with the reciprocal plasma-/? parameter 

Pmag 



Pr 



(16) 



(17) 



is known to be an important parameter characterizing the mor- 
phological appearance of jets with purely toroidal magnetic field 
configurations (where we follow the definitions from |Leismann 
|et al.|(|2005| )). For such jets, the pioneering work by Komissarov 
( |1999| l lias shown that Poynting flux dominated (high cr), 
strongly magnetized (high B r ) jets may develop a sharp 'nose 
cone', where highly magnetized plasma accumulates beyond the 
terminal shock of the jet beam. In contrast, strong backflows and 
therefore cocoon-dominated jets emerged when the kinetic en- 
ergy flux is dominant (low cr). In the |Leismann et aL ( 2005| l pa- 
rameter study extending these results to also purely poloidal field 
configurations, no 'nose cones' were found for poloidal topolo- 
gies. 

Due to the 2.5D nature of the equilibrium, all these param- 
eters in fact vary mildly to strongly across the beam radius. In 
Table [2] we computed beam-averaged values 



(18) 



nRf jo 



/=— 2 | InRfdR. 



From the values in this table, the sequence of models NR, MR, 
Ref 1 with increasing a clearly has similar inlet conditions (low 
cr, equipartition like fields with B r » 0.3), and explores the trend 
to f ss 7 flow regimes. In Figure [I] we plot the actual variation 
of the Lorentz factor F across the inlet for the Refl model, and it 
is seen that this model reaches axial flows with F up to 22. The 
figure also shows the inverse pitch p. for the same model. The 
inverse pitch fi is a quantitative measure of the twist in the jet 
magnetic field and is defined from 



~R~fh 



(19) 



This profile is important for a stability analysis of helically mag- 
netized jet flows, and as we have axial values of this inverse 
pitch of order 0.3 and as the current is distributed over the 
entire jet section ( |Appl et al. 1999]), we may expect that non- 
axisymmetric current-driven instabilities are less likely to play a 
role in the dynamics. However, a quantitative answer to stability 
with respect to non-axisymmetric perturbations (both of Kelvin- 
Helmholtz as well as current-driven type) is presently lacking 
for radially structured jets such as those introduced here. While 
simplified 'top-hat' profiles (i.e. layers of uniform, axially mag- 
netized flows) can still be analyzed analytically (Hard ee|2007[ l, a 
numerical analysis of the linearized relativistic MHD equations 
is needed to predict 3D effects which we artificially suppress 
by axisymmetry assumption. Naturally, the sequence of models 
Pol, Tor, Ref2 with varying field structure will be characterized 
by very different stability properties against non-axisymmetric 
modes. 

As a means to anticipate the 2.5D results, and to verify 
whether the employed resolution suffices to capture the flow de- 
tails, we first show the result of solving two Riemann problems 
that relate to the initial conditions in the reference case (Refl). 
We solve these Riemann problems numerically using essentially 




o.o 

0.0 0.5 1.0 

R 



Fig. 1. Jet inlet profiles for the reference model Refl . Left: radial 
variation of the Lorentz factor F. Right: variation of the inverse 
pitch p. 

Table 1. Model input parameters, as appearing in the profiles 
quantifying the radial variation of primitive variables. 



Model 


Pj/Pc 


Bo, B c , B\ 


a 


Pi 


NR 


0.1 


1,0.01, 1 


4.99 




MR 


0.1 


1,0.01, 1 


9.9 




Refl 


0.1 


1,0.01, 1 


9.99 




Pol 


0.1 


1,0.01,0.01 


999.9 




Tor 


0.1 


0, 0, 1 


9.99 




Ref2 


0.1 


2, 0.01, 1 


9.99 




Ref3 


0.2 


1,0.01, 1 


9.99 




Ref4 


0.1 


1,0.01, 1 


9.99 


2 



Table 2. Model parameters, in terms of beam-averaged inlet val- 
ues for various dimensionless ratios. 



Model 


r 


Pr 


cr 


M,M 


M a 


NR 


1.15 


0.3015 


0.0067 


3.6, 4.1 


6.6 


MR 


4.79 


0.2902 


0.0064 


7.1, 34.2 


13.6 


Refl 


6.92 


0.2899 


0.0064 


7.2, 50.3 


13.7 


Pol 


7.78 


0.2833 


0.0063 


7.2, 56.9 


13.9 


Tor 


6.92 


0.0009 


0.00002 


7.8, 53.3 


664.5 


Ref2 


6.92 


0.9147 


0.0250 


6.2, 44.7 


7.1 


Ref3 


6.92 


0.2899 


0.0064 


7.2, 50.3 


13.7 


Ref4 


6.92 


0.1514 


0.0064 


5.3, 36.5 


13.9 



identical settings for grid refinement and resolution. The on-axis 
conditions, where Lorentz factor r « 22 conditions prevail, are 
as in the ID problem where for constant B x = 1, left state LI is 
adjacent to right state Rl given by 





P 


P 


V x 


Vy 


By 


LI : 


100 


1 


0.999 








Rl : 


1000 


1.5 











L2 : 


100 


1 


0.992 


0.015 


0.15 


R2 : 


1000 


1.5 












(20) 



States L2 and R2 give a Riemann problem that relates to the 
conditions at radius R « 0.75, midway through the jet. It must 
be noted that this relation is somewhat ad hoc, since the two- 
dimensional case has significant variation of Bz(Z) and p(Z) at 
the head of the jet, which can not be mimicked in a ID Riemann 
problem. Nevertheless, Figs. 



3 show the result at t » 200 



for these Riemann problems. We capture in all cases both fast 




Fig. 2. A ID Riemann problem indicative for the on-axis condi- 
tions of the reference jet. Shown is the proper density, pressure, 
and Lorentz factor, at time t m 200. 



Fig. 3. A ID Riemann problem indicative for R * 0.75 (midway) 
conditions of the reference jet. Shown is the density, pressure, 
Lorentz factor, and v x , v y and B y at time t * 200. 



forward and reverse shocks very accurately, despite the fairly 
extreme parameters. The contact discontinuity is smeared out 
over many cells, but thanks to the AMR, adequately resolved. 
More problematic is the separation between contact discontinu- 
ity and the varying tangential field and velocity components in 
the Alfven signals in the second Riemann problem. Still, their 
amplitude and variation is appropriately represented (and the en- 
tropy is constant through these discontinuities). These ID tests 
qualify the shock-related pressure variations to be of several or- 
ders of magnitude, and indicate that we will be able to follow jet 
dynamics on a similar timescale within the computational do- 
main. 

4.2. Jet morphologies: from non-relativistic to relativistic 

To qualify the overall jet morphologies, we use the sequence of 
models characterized by increasing velocity (NR, MR, Refl). In 
the near non-relativistic model NR, the jet beam terminates in a 
Mach disk, across which a fair part of the directed kinetic energy 
is transferred to internal energy. This gives rise to a prominent 
hot spot of high pressure material, typically extending up to the 
contact discontinuity (or working surface) between shocked jet 
material, and shocked cloud material. The magnetic field is com- 
pressed by the reverse shock (Mach disk), enhancing the col- 
limation efficiency which limits the sideways expansion of the 
shocked beam matter (in between Mach disk and working sur- 



face). The shocked cloud material is bounded by an overarching 
bow shock. In this NR model, the inertia ratio between the jet 
and the external medium is low rj ~ f^pb/pc = 0.132, so that 
we expect a turbulent behavior near the working surface, which 
will lead to fairly strong disturbances to the jet. The compression 
of the shocked external medium is moderate, and the sideways 
spreading of shocked beam material is reflected by the external 
medium leading to a bow shock and backflow. In animated views 
of the jet evolution, one can witness the formation of vortical pat- 
terns at the contact interface in front of the Mach disk, which get 
pushed away laterally and form complex backflow patterns sur- 
rounding the jet beam proper. In doing so, recurrent cross-shocks 
are driven into the beam in the region behind the terminal shock. 
These lead to interacting diagonal cross-shock patterns internal 
to the beam, and they continually restructure the shape and lat- 
eral extent of the Mach disk at the end of the beam. A plot at 
time t =t 555 of the pressure distribution (showing the high pres- 
sure hot spot), the reciprocal plasma parameter /3 r , and the proper 
density (in logarithmic scale) is shown in Fig. [4] 

The magnetic field configuration for model NR at the same 
time t = 555 is shown in the top panel of Fig. [5] In the beam, 
slight helicity changes are associated with each diagonal cross- 
schock pattern, with increased twist regions across a converg- 
ing cross-shock, up to the diverging fronts. In the cocoon region 
formed by the backflowing material, rather strong azimuthal 
field components prevail, and the helically magnetized jet beam 
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Fig. 4. For the low speed T ^ 1.15 jet (NR), we show the pressure and reciprocal plasma beta j3 r (left) and proper density (right) 
distribution at time t =t 555, corresponding to 185 light crossing times of the jet diameter. 



gets surrounded by predominantly toroidal field regions. This 
can be seen in Fig. [5] where the right panel gives a view of se- 
lected field lines, combined with a translucent surface contain- 
ing a Schlieren plot of the rest frame density p. In the left panel, 
the velocity field shows the vortical motions associated with the 
backflows in the cocoon. The right half of this panel maps the 
distribution of the absolute value of the inverse pitch | fi | from 
Eq. 19 Regions where | fj. | exceeds unity are colored in black, 



and they give an indication of the most strongly wound field re- 
gions. It is seen how the p variation across the inlet is rather 
preserved throughout the jet beam up to the Mach disk, with 
modest variations associated with the diagonal cross-shocks as 
mentioned above. Beyond the Mach disk and in the cocoon, pre- 
dominant toroidal field prevails, Hence, the reverse shock (Mach 
disk) compresses the toroidal magnetic field. There is also clear 
evidence for a more poloidal, straight field layer in between the 
jet beam proper and the backflow regions. These backflows push 
the poloidal magnetic field towards the axis, and the resulting 
layer of strong poloidal magnetic field can be expected to in- 
crease the lateral stability of the jet beam. 

At all times in the evolution, only the jet beam, cocoon re- 
gions, and compressed region between beam and backflow, con- 
tain any significant magnetic pressure. This is shown in Fig. [6] 
where we now contrast the three jet models NR, MR, and Ref 1 
near the end of the simulated time intervals on the entire com- 
putational domain. Note that the non-relativistic model is shown 
at t ^ 555, while the faster jets MR and Refl are plotted at 
near identical, earlier times t ^ 220. The bottom half of each 
panel quantifies the magnetic pressure distribution p m . dg (in a 
logarithmic scale), and dynamically important values are found 
only up to the (turbulently deformed) contact interface between 
shocked beam plus backflow material with shocked cloud matter. 
This figure also shows that the faster jet models have a clearly 
more elongated appearance, consistent with their higher inertia 
oc T 2 (p + ypliy - 1)). The same figure also demonstrates that all 
our computations completely resolve the full bow shock pattern 
as they do not (yet) cross the lateral boundary. In figures through- 
out this document, we will show appropriately scaled, zoomed 
in regions of the domain only, to better illustrate flow details. In 



Fig. [5] the lower panel shows the magnetic field and flow infor- 
mation for the reference model Refl, in the same manner as the 
top panel for the non-relativistic case NR. Note that times, ver- 
tical extent and aspect ratio are quite different between the two 
cases, as evident from comparing them with Fig. [6] In the fast 
jet Refl (and also in the MR jet, which is not shown in Fig. |6j, 
we again see that the magnetic field in the jet beam roughly pre- 
serves the inlet helicity profile up to the Mach disk. The beam 
is bordered by a near vertical field 'sheet', and the vortical pat- 
terns marked by lower density structures contain the most thigh- 
tly wound fields. Unlike the non-relativistic model though, the 
backflowing vortices do not really form a rather extended co- 
coon surrounding the jet beam, but appear as more isolated nar- 
row protrusions into the jet cavity bounded by the bow shock. 
The locations of the dominant toroidal field regions are there- 
fore also rather localized. Thin strands of high field twist thus 
coincide with the lowest proper density spots which mark the 
centers of the vortical patterns. As will be discussed later on, the 
rotational flow patterns in these vortices are locally supersonic. 

In the preceding discussion, we contrasted flow morphology 
and dynamics for mildly relativistic to strongly relativistic cases 
(NR to Refl) all obtained in axisymmetric computations. All 
impulsively injected jet studies (where the jet gradually enters 
the simulation domain) tend to evolve from a 'ID' to a '2D' 
phase once eddies start to dominate the cocoon dynamics and 
internal jet beam shocks fully develop. For the models discussed 
here, the first internal cross-shock has already clearly formed at 
time t m 15. For the Refl case, multiple internal cross-shocks, 
and complex (ring-like) vortex sheddings are prominent beyond 
times t * 60. A transition to a '3D' phase where jet disrup- 
tion and/or significant non-axisymmetric deformation occurs is 
precluded in our setup. According to the discussion on 2D ver- 
sus 3D effects in Carvalho & O'Dea ( 2002) 1, jet disruption for 
classical HD jets may occur beyond distances given by SRjMj, 
where Mj indicates the jet Mach number. If we take the jet aver- 
aged Mach number M from Table [2] disruption to 3D can occur 
beyond time t » 54 and is likely to affect the further propaga- 
tion. Using the more appropriate relativistic value M ~ TM, 
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Fig. 5. For the NR jet from Fig.|4]at time t 555 on domain 20 x 60, the top panels show the flow topology on top of density plot 
(leftmost panel), along with a quantification of the inverse pitch throughout the jet: values above fi = 1 are colored black. In the right 
panel, a translucent Schlieren plot of the density is combined with field lines. The latter are colored according to their radial starting 
position (progressivley more outwards: blue R = 0.3, green R = 0.6, red R = 0.9, yellow R — 1.3, black R = 1 .8). The bottom panel 
gives the same info for the reference model Refl, at a different time t 220, using a very different aspect ratio of 20 x 160. 




Fig. 6. The (logarithm of) magnetic pressure p m . dg (bottom) and proper density p (top) for models NR (at t ^ 555), and MR and Ref 1 
(at t 220). The aspect ratio is identical in all frames, namely 40 x 200. 
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Fig. 7. The temporal evolution of the jet head position, for the sequence from non-relativistic, to the reference case (left panel), and 
a comparison with an approximate formula for kinetic energy dominated jets. At right, all simulated models are shown in a similar 
fashion, essentially demonstrating similar propagation speeds. 



this disruption length still lies beyond the computed time inter- 
val though. Earlier 3D non-relativistic hydrodynamic (e.g. Bodo 



(2002 



et al. (|1998|l) and mag netohydrodynamic (e.g. Baty & Keppens 
I; O'Neill et al. ( 2005) l) studies have been able to qualify 



some of the consequences of the axisymmetry restriction. Bodo 



|et al.| ( [T998[ ) identified how small-scale structure emerges fast in 
full 3D hydro simulations of periodic jet segments, and particu- 
larly light jets (just as those considered here) are rather prone to 
non-axisymmetric mode development. Dynamically important 
magnetic fields can help to stabilize the jet flow and mitigate the 
energy cascade to smaller scales. Classical 3D MHD simulations 



of periodic jet segments by Baty & Keppens (|2002) identified 
how helically magnetized jets can maintain jet coherency despite 
additional global non-axisymmetric instabilities. |Q'Neill et al. 
(2005 ) investigated propagation characteristics of helically mag- 
netized, light jets in 3D, following their propagation for lengths 
exceeding 100 jet radii. The presence of an intricately structured 
3D 'shock-web complex' at the frontal part of the jet was clearly 
seen in renderings of the compression rate, most prominent in 
the high Mach number jet propagating in uniform surroundings. 
In view of these results, it can be expected that full 3D relativistic 
simulations for the sequence NR to Ref 1 will show profound dif- 
ferences from our axisymmetric computations in the later stages 
of the simulation. 



4.3. Propagation speeds 

A straightforward qualification of the jet dynamics is obtained 
from their propagation speeds. It was shown in relativistic hy- 
dro ( Martf et al.|1997|l, and later used for relativistic MHD jets 
( |Leismann et al . 2005), that an estimate for the head advance 
speed could be made as 



v head - 



i + V^7I« 



vz 



(21) 



in which the jet beam internal versus (static F a = 1) ambient 
medium enthalpies ^/F 2 appear. This relativistic hydro estimate 
assumes pressure-matched conditions, and essentially uses mo- 
mentum balance in the shock frame. In our magnetized jets, 
all taken in the kinetic energy dominated regime (cr «: 1 see 
Table [2]), we can expect a reasonable agreement with this for- 
mula. Due to the internal beam equilibrium profiles, we can use 



this expression to get a 'lower' and 'upper' bound for the ex- 
pected jet propagation. Using averaged beam profile values in 
Eq. pT) , a reasonable lower bound is obtained, while when spec- 
ifying to beam axial values, a kind of higher bound is found. In 
Figure [7] the left panel shows the thus estimated speed ranges 
and compares them to the simulation results for the sequence 
of low to high speed jets (NR, MR, and Refl). For the slowest 
jet model, both estimates essentially coincide, and this jet shows 
good agreement with the predicted speed up to times t m 100, af- 
ter which the jet starts to propagate slower. For the faster MR to 
Refl models, there is also a trend to gradual deceleration, where 
the upper estimate prevails in the first phase of the computa- 
tion, but where we typically find actual propagation values in 
between the estimated bounds in the later stages. The right panel 
of Figure|7]compares the computed propagation velocities for all 
8 jets. The fastest jets in this limited sample are the nearly purely 
poloidal case, in accord with its higher beam (average as well as 
maximal, this time even going up to 70) Lorentz factor. In fact, 
this jet with (nearly) purely poloidal magnetic field behaves as 
a hydrodynamic jet (Majora ma & Anile|1987| l in first approxi- 
mation. Also the Ref3 case which has the reference jet structure 
penetrating a twice lighter external medium is fast, as it interacts 
weakly with the external medium. The propagation characteris- 
tics of the other four models turn out to be very similar, and this 
again agrees with their estimated values from Eq. pTl. 



Still, prominent 2D effects arise in our simulations, as e.g. 
seen in the shape of the bow shock in Fig. [6] This is influenced by 
the formation and interaction of the beam internal cross-shocks, 
as well as by the complex evolution of the final Mach disk. The 
first cross-shock for the NR case forms near Zsj4, and approx- 
imately remains at this location for times up to t w 300. For the 
MR case, the first cross-shocks initially forms at about Z ^ 12, 
gradually moving to larger distances. The reference Refl case 
initially forms a strong cross-shock at about Z « 25, but its lo- 
cation and strength varies with time: at t w 156, a fairly strong 
shock is visible at about Z * 64, but eventually a weaker shock 
remnant can be seen in Fig. [8] at about Z = 90 (followed by a 
much stronger shock at Z as 117). Despite these and other highly 
time-varying effects at the jet head, we find that the high Mach 
number jets considered here thus agree well with ID propagation 
estimates, consistent with findings for non-relativistic, axisym- 
metric hydro jets presented in |Carvalho & O'Dea| ( |2002"] l. 
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Fig. 8. For the reference Ref 1 case, we visualize the flow pattern in the cocoon, along with a local Mach number quantification 
(bottom half of top panel). Lower three panels: we show a cut along the symmetry axis of the Lorentz factor at the same time 
t =! 220, a cut of Br(Z) and fi^(Z) at the same time for R = 0.75. 



In Fig. [8] we illustrate several properties of the flow field for 
the Refl reference model. Using the same overall view as dis- 
played in Fig. [6] the bottom half of the top panel indicates where 
regions of locally supersonic poloidal velocities are encountered. 

All regions where ^jv^ 



+ v\jc s exceed unity are colored black, 



and the remaining greyscale is used for the subsonic range. The 
entire jet beam, together with the uppermost part of the bound- 
ing bow shock is locally supersonic. The bow shock speed varies 
significantly along its entire boundary. As its speed is mostly 
transversal, we indeed confirm that its lateral supersonic expan- 
sion varies strongly behind the jet head as already found from ax- 



isymmetric hydro studies ( Carvalho & O , Dea|2002| l. The upper 
half of this panel in Fig. [8 shows the (logarithm of) the Lorentz 
factor in greyscale, and it is seen that only the jet beam up to the 
Mach disk is seen to travel at a significant fraction of the speed 



of light (up to on axis values F w 22, with some variation across 
the internal cross-shocks). The velocity vectors are only drawn 
throughout the surrounding cavity, and show the complex vor- 
tical patterns, the overal expansion of the cavity, and clear evi- 
dence of localized regions of strong velocity shear. The latter aid 
in the formation of the turbulent structures by Kelvin-Helmholtz 
instabilities. The lower panels in Fig.|8]quantify the on-axis vari- 
ation of the Lorentz factor at this endtime and the magnetic 
field variation midway the jet beam (at radius R = 0.75). This 
shows that the beam internal cross-shocks act to convert some 
of the directed kinetic energy, and due to slight variations of the 
cross-section of the beam, matter accelerates again up to the next 
cross-shock. The magnetic field plays an important role in this 
process, as the field strongly pinches the flow downstream of the 
shock. Note how in this reference case, Lorentz factor T ^ 20 
flow occurs up to distance of Z ^ 158, while at the final Mach 
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Fig. 9. For the reference jet Refl, we show cuts through the jet head, taken at fixed radius R — 1 and for Z e [161, 171], and at 
the time t m 220. From left to right and top to bottom, we show field components B R , B z , B v , proper density p, Lorentz factor F, 
pressure p, along with entropy S , and velocity components vz and v v . 



disk, the Lorentz factor is still above 10. In order to quantify 
the deceleration process better, even longer simulations would 
be required. 

5. Jet head structure and varying magnetic 
topology 

If we turn to the jet head structure in particular, we can quan- 
tify the compression occuring at the Mach disk, and confront 
the field structure with that obtained in the ID Riemann prob- 
lem shown in Fig. [3] Fig. [9] shows the axial variation at a radial 
distance R = 1, zoomed in on the jet head structure (note the lim- 
ited Z-range). We discuss this structure from front to back. Then, 
one first detects the forward shock (bow shock, at Z m 170.8) 
which leads to pronounced increases in pressure and density, of 
a magnitude consistent with what the ID problem demonstrated. 
Another discontinuity of a pure thermodynamic nature trails at 
Z * 169.9. In a zone Z e [169.6, 169.9], we find the highest 
values for proper density p m 6000 and pressure p * 1000. 
Magnetic field variations all occur prior to Z w 169.6, where 
the density and pressure decrease from the values mentioned, to 
values still significantly above their beam values. We can detect 
distinct variation of all three magnetic field components within 
the region Z e [168.8, 169.6], coinciding with increased entropy 
S = pp~ y within this zone. At the trailing shock (Mach disk) 
located at Z m 168.8, the Lorentz factor drops from about 6 to 2. 
The discontinuity seen at Z « 162 is the last internal shock prior 



to the jet head. There, as well as at the Mach disk, the azimuthal 
field is enhanced, in accord with the increase in transverse mag- 
netic field seen in the ID problem from Fig. [3] The analogy can 
not be expected to hold beyond such qualitative features: the 
actual 2D jet shows significant multi-D variation: the magnetic 
field becomes nearly purely radial and azimuthal within a thin 
layer Z € [169.3, 169.9], and for Z > 169.9 the uniform vertical 
cloud field of magnitude B c = 0.01 is found. 

The sequence of models termed Pol, Refl, Ref2, and Tor 
can be interpreted as a series with similar inlet characteristics in 
terms of average Lorentz factor and Mach number, and changing 
magnetic topology from near poloidal to pure toroidal. In accord 
with the decreasing influence of poloidal field in the jet beam 
(from Pol to Tor), the backflows tend to show more fine structure 
with more eddies being shed in the purely toroidal case. There is 
then accordingly more internal beam structure: multiple cross- 
shocks are induced and interact. Once more, full 3D simulations 
are needed to assess the generality of this result: the ringed vor- 
tex structures we obtain may break down and get stretched in 
non-axisymmetric fashion. This can alter the overall field topol- 
ogy drastically and impact the backflow. We find that in all but 
the purely toroidal case, the regions of significant magnetic pres- 
sure are in the beam and backflows. The purely toroidal case 
carried a weak magnetic field, and this can be seen to trace the 
region occupied by shocked beam matter forming the inner co- 
coon. 
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Fig. 10. For a sequence of models from toroidal (top), to a model with double the poloidal field magnitude as the reference case 
(middle), to nearly purely poloidal magnetic field (bottom): maps of arbitrary scaled 'power' given by T 2 v 2 B 2 sin 2 ( v F) (left), together 
with logarithmic greyscale plots of the temperature and pressure distributions. All plots correspond to instantaneous values at the 
endtime of the simulation, and only show a zoom on the top half of the entire domain. 



As a qualitative means to address how different magnetic 
topologies lead to varying jet morphologies we use this model 
sequence as follows. A formula quantifying the total emitted ra- 
diation by a single electron traveling at speed v in its relativistic 
cyclotron motion about a magnetic field B is given by (Rybicki 
& Lightman|T979] l 

2e 4 

Put = ^v 2 r 2 B 2 sin 2 ^) (Js- 1 ) , (22) 



\2neam 2 c 



where *F indicates the pitch angle between particle velocity v and 
magnetic field B (in Tesla). In our relativistic MHD simulations, 



we treat the plasma as a single fluid characterized by its bulk 
speed v. A crude means to infer observational 'synchrotron' in- 
tensity uses the bulk plasma speed in Eq. ( 22 1 together with the 



detailed knowledge on the magnetic field distribution. In Fig. 10 
the left panels produce arbitrary scaled maps of the essential lo 
cal dependence v 2 Y 2 B 2 sin 2 ^) for the jet sequence, at the end 
time of the simulations. It is to be noted that we present instanta 
neous local values of this 'power map', which is not appropriate 
for quantifying true synchrotron emission for these relativistic 
jets, as we e.g. do not incorporate time retardation effects. A 
procedure to incorporate all relevant relativistic effects uses the 
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Fig. 11. For a heavier jet (top), versus a higher pressure jet (bottom): the proper density distribution at the endtime of the simulations. 
A reduced density contrast ensures fairly stable jet flows and little internal as well as cocoon structure. Increased jet pressure leads 
to more internal beam structure. 



output of relativistic HD simulations and solves the relativistic 
radiative transfer equation assuming optically thin conditions, 
as done in [Komissarov & Falle| ( |1997| >. Recently, Zakam ska et| 



al. 



2008) used analytic self-similar axisymmetric RMHD jet 



models, to quantify both synchrotron emission as well as po- 
larization maps by performing the proper integration along the 
line of sight in the observer's reference frame. In this paper, we 
restrict ourselves to quantifying the 'power' dependence from 
Eq. p2) , since for practical purposes one would need to perform 
the data processing during the simulation, and make reasonable 
assumptions about the relativistic particle distribution function 
(which is not contained in the ideal RMHD simulation). The 
local instantaneous values give quantitative information on the 
overall flow and field topology, and it is this aspect which con- 
cerns us here mostly. Only about the front half of the simulated 
domain is shown in Fig. 10 and a simultaneous quantification 
of the pressure and 'temperature' T = pip is shown at right. 



The highest pressure and temperature regions coincide with the 
(very narrow) regions between the Mach disk and the contact 
surface, and show up as bright regions ahead of the jet beams. 
Consistent with the higher pressure found in the frontal com- 
pressed beam matter when more poloidal field is present, the 
associated (narrow) hot spot eventually dominates. In all cases 
though, the complex shock interactions seen in the beam clearly 
show up. Given the aforementioned trend towards more complex 
shock interactions when more toroidal field structure is present, 
the power maps reflect this trend directly. The angle-dependent 
factor in formula ( 22 1, together with the magnetic field variation 
with radius, also explains the trend seen from more outer beam 
sensitivity in the almost toroidal case, to pronounced inner beam 
sensitivity at diagonal cross-shock fronts for more poloidal field 
configurations. This hints at the possibility to infer field topol- 
ogy characteristics from true synchrotron emission maps, which 
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need to be constructed rigorously from simulations such as these 
presented here. 



6. Thermodynamic variations 



The final two models are shown in Figure 1 1 showing the loga 



rithm of their proper densities. The top panel is for the case with 
reduced density contrast, where the jet penetrates less dense ma- 
terial than in the reference case Refl. The magnetized jet travels 
correspondingly faster, and a remarkably stable beam with little 
backflow features results. The stand-off distance to the first in- 
ternal beam cross-shock is also increased. In contrast, the Ref4 
model which has a higher overall pressure compared to the ref- 
erence case, shows a similarly rich pattern of vortex structures. 
Since the magnetic field is identical in these models, the jet ap- 
pearance seems to relate most dramatically on density differ- 
ences, at least for these near-equipartition, kinetically dominated 
jets. Comparing the increased pressure jet case with the refer- 
ence case in Figure[6] the richly structured backfiows are slightly 
more pronounced, and the internal beam density contrasts are 
more evident. Within the simulated time, the higher pressure jet 
propagation has had about one more phase of slight decelera- 
tion, reformation of the Mach disk, and consecutive acceleration 
shaping the overall bowshock. This is in accord with its shorter 
internal sound crossing time. Quantitatively speaking, the high- 
est proper density values (always found in the front of the jet 
head) are found for the reference case Refl where its instanta- 
neous maximum is found p max = 6280 (see Fig. [6]), while the 
higher pressure case demonstrates values up to p max = 5236. For 
the lower density contrast, we find, as expected due to the lesser 
contrast and cloud density, a reduced value p max = 2818. 



7. Outlook 

We studied the morphology and propagation characteristics of a 
series of highly relativistic, helically magnetized jets. All were 
kinetic energy dominated, and apart from the series of increas- 
ing averaged beam velocity, were similar in propagation speeds 
and overall magnetization (equipartition). With modest changes 
in magnetic topology, as well as in internal pressure and den- 
sity, fairly distinct differences can be detected in the distribution 
of local power, as in the cocoon and internal beam structure. 
The radially stratified jets studied here show significant varia- 
tion of Lorentz factor across their diameter, and the high cen- 
tral 'spine' (exceeding F = 20) results in fairly elongated bow 
shocks. The magnetic helicity changes at internal cross-shocks 
act to reaccelerate the jet repeatedly. The lighter jets show more 
fine structure in their backfiows, and this will likely continue 
over larger density contrasts than those studied here. Previous 
studies (Leismann et al. 2005; Komissarov| |1999| ) considered 
density ratios of 100 or more. The models here only considered 
density ratios of 5 to 10, and investigate jet propagation in less 
dense environments which are known to lead to more stable jet 
configurations over longer distances. The sheet of more poloidal 
field surrounding the jet beam proper which we found in the sim- 
ulations also aids in stabilizing the beam. 

In our simulations, the resulting jet propagation exhibits a 
small deceleration of the head of the (still) relativistic jet along 
its axis, only seen in the decrease in Lorentz factor along the axis 
occuring at the internal shocks. Such deceleration has been ob- 
served in various FR II jets, in correlation with an increase of the 
magnetic field intensity and particle density ( Georgeanopoulos 
|& Kazanasl [20041 |Sambruna et~aT1 |2006b| >. |Tavecchio et~aT 



(2006) have shown that this deceleration is compatible with the 
entrainment of the external gas by the jet. We plan to perform 
longer term simulations, in order to obtain a better description 
of the jet braking as well as a study of the impact of the external 
mass density compared to the jet one. In this kind of computa- 
tion, we will explore regimes where the jet kinetic energy and 
the jet electromagnetic energy are of the same order. The use of 
the AMR strategy will be a powerful tool in the study of the as- 
sociated magnetic field and density amplification occurring near 
the jet head. We then aim to provide statements from modeling 
results, about the ratio of the jet density to the external medium 
density, and use it to constrain better the observational values. 

We also intend to explore the computationally challenging 
regime of Poynting flux dominated jets with helical field topolo- 
gies in future work. Since most AGN jets are likely associated 
with high Poynting flux jets near the source, it will be of interest 
to study the transition from Poynting to kinetic energy domi- 
nated jets in even larger scale computations. Also, the stability 
of these jets in fully three-dimensional simulations is as yet un- 
explored. We demonstrated that the inlet jet magnetic topology 
for the sample studied here seems to be maintained over a sig- 
nificant distance, and the helicity measure showed only strong 
toroidal field concentrations in the localized vortices developing 
from the backfiows. This was consistent with a trend followed 
from slower 'non-relativistic' jets, where toroidal field gets cre- 
ated across the Mach disk, and subsequently mixed into the co- 
coon. Since the toroidal field concentrations for the reference 
helical jet are mostly in localized, supersonically rotating vor- 
tices, their tendency to induce kink deformations may be less. 
Extremely high resolution (grid-adaptive) 3D studies are called 
for to investigate this issue further. 
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